/* STANDARD LIBRARY */
#include<cstdio>
#include<stdlib.h>
#include<iostream>
#include<string>
#include<math.h>
#include<ctime>

/* ADDITIONAL LIBRARY */
#include"constants.h" // constants
#include"supToGSL.h"

/* DECLATATIONS */
#ifndef DECLARATIONS_supToJK
#define DECLARATIONS_supToJK

/* Gutzwiller factor from Gutzwiller paper
 * my own calculation usinig Mathematica */
class Gutzwiller
{
public:
	// for definitions of functions see supToMain.cpp file
	double value(double n, double m, double d, double sigma);
	double dn(double n, double m, double d, double sigma);
	double dm(double n, double m, double d, double sigma);
	double dd(double n, double m, double d, double sigma);
};

/* Dispersion relation*/
class Dispersion
{
protected:
	const double t1; // hopping integral, the nearest neighbors
	const double t2; // hopping integ., second the nearest neighbors
public:
	Dispersion(const double _t1, const double _t2) : t1(_t1), t2(_t2) {}
	double value2D(double kx, double ky)
	{
		double cosx = cos(kx);
		double cosy = cos(ky);
		return -2.*t1*(cosx + cosy) + 4.*t2*cosx*cosy;
	}
	virtual double value3D(double kx, double ky, double kz)
	{
		double cosx = cos(kx);
		double cosy = cos(ky);
		double cosz = cos(kz);
		return -2.*t1*(cosx+cosy+cosz) + 4.*t2*(cosx*cosy+cosx*cosz+cosy*cosz);
	}
};

class DispersionZdependence : public Dispersion
{
protected:
	const double t1z;
	const double t2z;

public:
	DispersionZdependence(const double _t1, const double _t2, const double _t1z, const double _t2z) \
	: Dispersion(_t1, _t2), t1z(_t1z), t2z(_t2z) {}

	double value3D(double kx, double ky, double kz)
	{
		//fprintf(stderr, "DispersionZdependence: %f, %f, %f", kx, ky, kz);
		double cosx = cos(kx);
		double cosy = cos(ky);
		double cosz = cos(kz);
		return -2.*t1*(cosx + cosy) + 4.*t2*cosx*cosy - 2.*t1z*cosz + 4.*t2z*(cosx*cosz+cosy*cosz);
	}
};

/* Feermi Distribution */
class Distribution
{
public:
	double value(double x, double beta)
	{ return 0.5 * (1.0 - tanh( 0.5 * beta * x )); }
};

/* CLASSES */
class JK
{
protected:
	Gutzwiller *g;
	Dispersion *e;
	Distribution *f;
	double h, n, t1, t2, t1z, t2z, beta, U;
	unsigned int L;
	const unsigned int N;
	enum Variables { mu, m, d, l_n, l_m };
	double *y;

	void clear(double *t)
	{ for(unsigned int i=0; i<N; i++) { t[i]=0; } }

public:
	JK(Parameters *Par, int numberOfEquation);
	~JK();

	virtual void onePoint(double *x, double *y, int i, int j) = 0;
	virtual void onePoint3D(double *x, double *y, int i, int j, int k) = 0;
	virtual void OneLayer (double *x, double *y);
	virtual void OneCubic (double *x, double *y);
	virtual double* sumation2D(double *x) =0;
	virtual double* sumation3D(double *x) =0;
};

class JKGA : public JK
{
protected:
	enum Variables { mu, m, d };

public:
	JKGA(Parameters *Par) : JK(Par, 3) { }
	~JKGA() {}

	void onePoint(double *x, double *y, int i, int j);
	void onePoint3D(double *x, double *y, int i, int j, int k);
	double* sumation2D(double *x);
	double* sumation3D(double *x);
};

class JKSGA : public JK
{
private:
	enum Variables { mu, m, d, l_n, l_m };

public:
	JKSGA(Parameters *Par) : JK(Par, 5) { }
	~JKSGA() {}

	void onePoint(double *x, double *y, int i, int j);
	void onePoint3D(double *x, double *y, int i, int j, int k);
	double* sumation2D(double *x);
	double* sumation3D(double *x);
};

// for threads
class ParThre
{
public:
	int id; // id of current thread
	double *x0; // start point: // mu, m, d
	double *par; // parameters: h, n, t1, t2, beta, U, L
	unsigned int d; // which parameter will be changed
	unsigned int nRepetitions; // and how many times
	double start; // from which value
	double end; // to which value
	size_t nEquations; // how many equations
	unsigned int MaxNumberOfShake; // how many we can shake system?
	const char *file; // to which file write results?
	const char *title; // briefly describe - what is it?
	unsigned maxInterp; // how many points will be use to interpolation
	int (*func)(const gsl_vector*, void*, gsl_vector*);
	pthread_mutex_t *mut; // mutex
	unsigned int nPar; // how many parameters we have

public:
	ParThre() {}
};

#endif
